! Bessel Y_0(x) function in double precision
!
      function dbesy0(x)
      implicit real*8 (a - h, o - z)
      dimension a(0 : 15), b(0 : 111), c(0 : 125), d(0 : 51)
      parameter (pi4 = 0.78539816339744830962d0)
      data (a(i), i = 0, 15) / 
     &    -0.00000000000151249795d0, 0.00000000029979612902d0, 
     &    -0.00000004317352912436d0, 0.00000431735413787068d0, 
     &    -0.00027631066508933090d0, 0.00994718394324338940d0, 
     &    -0.15915494309189533339d0, 0.63661977236758134306d0, 
     &    0.00000000000409490035d0, -0.00000000076925095943d0, 
     &    0.00000010358472550303d0, -0.00000949500519343105d0, 
     &    0.00053860266685948738d0, -0.01607396802593822992d0, 
     &    0.17760601686906713536d0, -0.07380429510868722524d0 / 
      data (b(i), i = 0, 13) / 
     &    -0.00000000002958527319d0, -0.00000000004799424787d0, 
     &    0.00000000070761365903d0, 0.00000000786916310954d0, 
     &    0.00000004406321746330d0, 0.00000010474207054310d0, 
     &    -0.00000070359581447993d0, -0.00001068225180236166d0, 
     &    -0.00007636033201484949d0, -0.00034318294473105478d0, 
     &    -0.00050019381740765567d0, 0.00898034680285547482d0, 
     &    0.14783488851798565262d0, 0.01218096458136948754d0 / 
      data (b(i), i = 14, 27) / 
     &    -0.00000000003790394288d0, -0.00000000083012324083d0, 
     &    -0.00000000416344143065d0, -0.00000000159750317942d0, 
     &    0.00000011897479804282d0, 0.00000089462034665663d0, 
     &    0.00000279277855080310d0, -0.00000710267818373816d0, 
     &    -0.00014158631248747538d0, -0.00089512296835766028d0, 
     &    -0.00286419212921734238d0, 0.00448696141199596815d0, 
     &    0.16247276853497028694d0, 0.16806523426268420517d0 / 
      data (b(i), i = 28, 41) / 
     &    0.00000000049508820408d0, 0.00000000240256475451d0, 
     &    -0.00000000341287024760d0, -0.00000007926555593644d0, 
     &    -0.00000027510687602104d0, 0.00000082850014912691d0, 
     &    0.00001202838962100801d0, 0.00004446440998573947d0, 
     &    -0.00006345799532264893d0, -0.00153636349657759957d0, 
     &    -0.00784582749467635199d0, -0.01091113270659651752d0, 
     &    0.15855078058220179592d0, 0.33112075987570150227d0 / 
      data (b(i), i = 42, 55) / 
     &    -0.00000000121381163034d0, -0.00000000185252445637d0, 
     &    0.00000002782181715752d0, 0.00000009983716245339d0, 
     &    -0.00000051382328105461d0, -0.00000413045112612025d0, 
     &    0.00000069219807212556d0, 0.00011346881914173059d0, 
     &    0.00044857142622311310d0, -0.00075862241142011713d0, 
     &    -0.01330170710405210993d0, -0.04337207915110066583d0, 
     &    0.10708263075506099964d0, 0.46937172516119821029d0 / 
      data (b(i), i = 56, 69) / 
     &    0.00000000123954313286d0, -0.00000000131506473470d0, 
     &    -0.00000003208158247443d0, 0.00000006936208476862d0, 
     &    0.00000096326321605564d0, -0.00000186813677726279d0, 
     &    -0.00003296596487362508d0, -0.00001027532398838717d0, 
     &    0.00088296805857821394d0, 0.00288502916411875494d0, 
     &    -0.00981551522729451135d0, -0.08175499367255420548d0, 
     &    -0.01970902432425660780d0, 0.51958016003260732459d0 / 
      data (b(i), i = 70, 83) / 
     &    -0.00000000050293922567d0, 0.00000000259058191150d0, 
     &    0.00000000380069015768d0, -0.00000012579287605830d0, 
     &    0.00000029666593806065d0, 0.00000532199410651341d0, 
     &    -0.00001466013674557952d0, -0.00021196179800800308d0, 
     &    0.00014889075157975250d0, 0.00598670133725335804d0, 
     &    0.00917434640154761191d0, -0.08592436550224288681d0, 
     &    -0.19700894914535846940d0, 0.41202451505149651586d0 / 
      data (b(i), i = 84, 97) / 
     &    0.00000000002689416111d0, -0.00000000060225798842d0, 
     &    0.00000000716632287642d0, -0.00000003038834206093d0, 
     &    -0.00000049829955736697d0, 0.00000364678661079850d0, 
     &    0.00002599276311771812d0, -0.00016374249229498868d0, 
     &    -0.00112385259950125883d0, 0.00342431039748960465d0, 
     &    0.03016663018029814938d0, -0.02432512430353342927d0, 
     &    -0.31797371916576712163d0, 0.14418001571733803493d0 / 
      data (b(i), i = 98, 111) / 
     &    0.00000000000194632210d0, -0.00000000013451327719d0, 
     &    0.00000000336581427068d0, 0.00000002163920016772d0, 
     &    -0.00000050868502638387d0, -0.00000126562372408989d0, 
     &    0.00003554054033687039d0, 0.00007490048474688881d0, 
     &    -0.00142380042008861330d0, -0.00355314274445060456d0, 
     &    0.03042001885994619102d0, 0.07365497405664882878d0, 
     &    -0.26882214651298246524d0, -0.16578636480856581212d0 / 
      data (c(i), i = 0, 13) / 
     &    0.00000000014133441141d0, -0.00000000116239182574d0, 
     &    0.00000000560813675548d0, 0.00000004361111919906d0, 
     &    -0.00000019980743680592d0, -0.00000623828267379045d0, 
     &    0.00002659779111113435d0, 0.00034429768949740887d0, 
     &    -0.00128058451746064875d0, -0.01165906901727504569d0, 
     &    0.03800023545011044967d0, 0.13079665132249175616d0, 
     &    -0.30099732306965462304d0, -0.19470500862950453349d0 / 
      data (c(i), i = 14, 27) / 
     &    0.00000000002985070388d0, -0.00000000043935370028d0, 
     &    -0.00000000270094988368d0, 0.00000005270139367180d0, 
     &    0.00000035596645009533d0, -0.00000559393055101660d0, 
     &    -0.00002415386269153308d0, 0.00034989509210648774d0, 
     &    0.00098204687587121096d0, -0.01241793885175385487d0, 
     &    -0.01398519840974389960d0, 0.16758045653582919005d0, 
     &    0.02375823895638961834d0, -0.33948059288191103829d0 / 
      data (c(i), i = 28, 41) / 
     &    0.00000000003440543098d0, -0.00000000001492311224d0, 
     &    -0.00000000550307983816d0, 0.00000000288861169685d0, 
     &    0.00000065975576340131d0, -0.00000064253433004282d0, 
     &    -0.00005095354686815779d0, 0.00006355454993858760d0, 
     &    0.00231750171436134067d0, -0.00344140664586595266d0, 
     &    -0.04796153689875432703d0, 0.06553727330877767204d0, 
     &    0.27409127395927545297d0, -0.17324243491898233567d0 / 
      data (c(i), i = 42, 55) / 
     &    0.00000000001869953108d0, 0.00000000035678205287d0, 
     &    -0.00000000323369443017d0, -0.00000004932837601739d0, 
     &    0.00000040636932259569d0, 0.00000455294558480472d0, 
     &    -0.00003376596464782138d0, -0.00025758322585307069d0, 
     &    0.00167416927479443727d0, 0.00735300543429220468d0, 
     &    -0.03904554680817849396d0, -0.07593187710651205994d0, 
     &    0.25912851048611625179d0, 0.11731328614820863082d0 / 
      data (c(i), i = 56, 69) / 
     &    -0.00000000001168499781d0, 0.00000000040491570620d0, 
     &    0.00000000176042185153d0, -0.00000005796810963324d0, 
     &    -0.00000017664264701199d0, 0.00000565232278687159d0, 
     &    0.00001060939366068237d0, -0.00034382320897779425d0, 
     &    -0.00028786871040915476d0, 0.01103696602501040460d0, 
     &    0.00105742483935856071d0, -0.13664188676516064192d0, 
     &    0.02616867939853747003d0, 0.27020510536578747599d0 / 
      data (c(i), i = 70, 83) / 
     &    -0.00000000003097140650d0, 0.00000000010275535029d0, 
     &    0.00000000507546062007d0, -0.00000001697584561546d0, 
     &    -0.00000058222618994516d0, 0.00000192483899143098d0, 
     &    0.00004338881805420200d0, -0.00013714074807064326d0, 
     &    -0.00184722933402281351d0, 0.00517360719982082990d0, 
     &    0.03611657815299529568d0, -0.07491163418624572802d0, 
     &    -0.20317989938720766824d0, 0.17121062620272384486d0 / 
      data (c(i), i = 84, 97) / 
     &    -0.00000000002289344046d0, -0.00000000027762757821d0, 
     &    0.00000000391330214230d0, 0.00000003673304032817d0, 
     &    -0.00000047258696909223d0, -0.00000323109113206792d0, 
     &    0.00003749284946160478d0, 0.00017039067180362504d0, 
     &    -0.00172637004639815403d0, -0.00454152531086626043d0, 
     &    0.03717220530377418124d0, 0.04489395902785574534d0, 
     &    -0.23370422835726857839d0, -0.06753037249787639679d0 / 
      data (c(i), i = 98, 111) / 
     &    0.00000000000473299038d0, -0.00000000040489249115d0, 
     &    -0.00000000056747261733d0, 0.00000005656135976496d0, 
     &    0.00000003564150059329d0, -0.00000534821220288168d0, 
     &    0.00000007581385461348d0, 0.00031190256463318867d0, 
     &    -0.00014634564376714538d0, -0.00958213551550479850d0, 
     &    0.00624681472076521329d0, 0.11513529702572439703d0, 
     &    -0.05794254714300082167d0, -0.22523211169118786537d0 / 
      data (c(i), i = 112, 125) / 
     &    0.00000000002732444065d0, -0.00000000017726475495d0, 
     &    -0.00000000437562701999d0, 0.00000002681536885560d0, 
     &    0.00000048799100503315d0, -0.00000276490187540644d0, 
     &    -0.00003513357925824063d0, 0.00017693155138571443d0, 
     &    0.00144533233498513085d0, -0.00599106092630544975d0, 
     &    -0.02759437856689911694d0, 0.07945362316083459396d0, 
     &    0.15383825653750118008d0, -0.17121430684466928734d0 / 
      data (d(i), i = 0, 12) / 
     &    1.059601355592185731d-14, -2.71150591218550377d-13, 
     &    8.6514809056201638d-12, -4.6264028554286627d-10, 
     &    5.0815403835647104d-8, -1.76722552048141208d-5, 
     &    0.16286750396763997378d0, 2.949651820598278873d-13, 
     &    -8.818215611676125741d-12, 3.571119876162253451d-10, 
     &    -2.631924120993717060d-8, 4.709502795656698909d-6, 
     &    -5.208333333333283282d-3 / 
      data (d(i), i = 13, 25) / 
     &    7.18344107717531977d-15, -2.51623725588410308d-13, 
     &    8.6017784918920604d-12, -4.6256876614290359d-10, 
     &    5.0815343220437937d-8, -1.76722551764941970d-5, 
     &    0.16286750396763433767d0, 2.2327570859680094777d-13, 
     &    -8.464594853517051292d-12, 3.563766464349055183d-10, 
     &    -2.631843986737892965d-8, 4.709502342288659410d-6, 
     &    -5.2083333332278466225d-3 / 
      data (d(i), i = 26, 38) / 
     &    5.15413392842889366d-15, -2.27740238380640162d-13, 
     &    8.4827767197609014d-12, -4.6224753682737618d-10, 
     &    5.0814848128929134d-8, -1.76722547638767480d-5, 
     &    0.16286750396748926663d0, 1.7316195320192170887d-13, 
     &    -7.971122772293919646d-12, 3.544039469911895749d-10, 
     &    -2.631443902081701081d-8, 4.709498228695400603d-6, 
     &    -5.2083333315143653610d-3 / 
      data (d(i), i = 39, 51) / 
     &    3.84653681453798517d-15, -2.04464520778789011d-13, 
     &    8.3089298605177838d-12, -4.6155016158412096d-10, 
     &    5.0813263696466650d-8, -1.76722528311426167d-5, 
     &    0.16286750396650065930d0, 1.3797879972460878797d-13, 
     &    -7.448089381011684812d-12, 3.512733797106959780d-10, 
     &    -2.630500895563592722d-8, 4.709483934775839193d-6, 
     &    -5.2083333227940760113d-3 / 
      if (x .lt. 0.85d0) then
          t = x * x
          y = ((((((a(0) * t + a(1)) * t + 
     &        a(2)) * t + a(3)) * t + a(4)) * t + 
     &        a(5)) * t + a(6)) * t + a(7)
          y = ((((((a(8) * t + a(9)) * t + 
     &        a(10)) * t + a(11)) * t + a(12)) * t + 
     &        a(13)) * t + a(14)) * t + a(15) + y * log(x)
      else if (x .lt. 4.5d0) then
          t = x - 4 / x
          k = int(t + 4)
          t = t - (k - 3.5d0)
          k = k * 14
          y = ((((((((((((b(k) * t + b(k + 1)) * t + 
     &        b(k + 2)) * t + b(k + 3)) * t + b(k + 4)) * t + 
     &        b(k + 5)) * t + b(k + 6)) * t + b(k + 7)) * t + 
     &        b(k + 8)) * t + b(k + 9)) * t + b(k + 10)) * t + 
     &        b(k + 11)) * t + b(k + 12)) * t + b(k + 13)
      else if (x .lt. 12.5d0) then
          k = int(x)
          t = x - (k + 0.5d0)
          k = 14 * (k - 4)
          y = ((((((((((((c(k) * t + c(k + 1)) * t + 
     &        c(k + 2)) * t + c(k + 3)) * t + c(k + 4)) * t + 
     &        c(k + 5)) * t + c(k + 6)) * t + c(k + 7)) * t + 
     &        c(k + 8)) * t + c(k + 9)) * t + c(k + 10)) * t + 
     &        c(k + 11)) * t + c(k + 12)) * t + c(k + 13)
      else
          v = 24 / x
          t = v * v
          k = 13 * (int(t))
          y = ((((((d(k) * t + d(k + 1)) * t + 
     &        d(k + 2)) * t + d(k + 3)) * t + d(k + 4)) * t + 
     &        d(k + 5)) * t + d(k + 6)) * sqrt(v)
          theta = (((((d(k + 7) * t + d(k + 8)) * t + 
     &        d(k + 9)) * t + d(k + 10)) * t + d(k + 11)) * t + 
     &        d(k + 12)) * v - pi4
          y = y * sin(x + theta)
      end if
      dbesy0 = y
      end
!
